! This sets the values of most of the key model parameters
module params

implicit none

integer,parameter :: prec=selected_real_kind(15,300)

integer ni, nt, nw, nwp, na, nm, nh, nk, nbignodes, nonodesk, nosims

!For baseline calibrated run, these "model structure" parameters are:
parameter (ni=2, nt=15, nw=6, nwp=3, na=20, nm=11, nh=3, nk=9, nosims=9999) 
!parameter (ni=3, nt=15, nw=6, nwp=3, na=20, nm=11, nh=3, nk=9, nosims=9999)  !for public option
!parameter (ni=2, nt=30, nw=6, nwp=3, na=20, nm=11, nh=3, nk=9, nosims=9999)  !for annual periods - also need to change ntHRS, yrsinpd, annualperiods, r, beta
!parameter (ni=2, nt=15, nw=6, nwp=3, na=20, nm=15, nh=3, nk=9, nosims=9999)  !for joint pareto-altruism

integer ierr,status,pid,nproc,node1,node2
integer whichiter
integer seed

integer, parameter :: DoEstimation=0
integer, parameter :: DoStandardErrors=0
integer, parameter :: DoIdentification=0

integer, parameter :: simlowfloor=0
integer, parameter :: simnoinf=0
integer, parameter :: simcashins=0
integer, parameter :: simcashmcd=0
integer, parameter :: simfullfairins=0
integer, parameter :: simpublicoption=0
integer, parameter :: noncoop=0
integer, parameter :: estimatesimple=0
integer, parameter :: annualperiods=0

integer, parameter :: Nparest=9
integer, parameter :: Nmoments=94
integer, parameter :: limcom=1
integer, parameter :: parallel=1

integer, parameter :: nachoice=30
double precision :: achoice(nachoice)

integer momi, ci
real(prec) :: estmatrix(10000,Nparest+1)
real(prec) :: moment_data_vector(Nmoments),moment_simu_vector(Nmoments),moment_err_vector(Nmoments)
real(prec) :: weight_matrix(Nmoments,Nmoments)
real(prec) :: varcov_matrix(Nmoments,Nmoments)

real(prec) formalpref
real(prec) mcaidpref
real(prec) stigma
real(prec) alpha
real(prec) gamma
real(prec) gammab
real(prec) termphiL
real(prec) termphiH
real(prec) altruism
integer    musimstartparam
integer, parameter :: musimstartparamc(4) =(/5,5,5,4/) !for sensitivity to starting cohorts at mean age-specific weight

real(prec) cminp
real(prec) cmink
real(prec) lshifter

real(prec) :: passet_c1_p25_data(9,2), kasset_c1_p25_data(9,2) !moment,wgt (inverse variance)
real(prec) :: passet_c1_p50_data(9,2), kasset_c1_p50_data(9,2) !moment,wgt (inverse variance)
real(prec) :: passet_c1_p75_data(9,2), kasset_c1_p75_data(9,2) !moment,wgt (inverse variance)
real(prec) :: passet_c2_p25_data(9,2), kasset_c2_p25_data(9,2) !moment,wgt (inverse variance)
real(prec) :: passet_c2_p50_data(9,2), kasset_c2_p50_data(9,2) !moment,wgt (inverse variance)
real(prec) :: passet_c2_p75_data(9,2), kasset_c2_p75_data(9,2) !moment,wgt (inverse variance)
real(prec) :: passet_c3_p25_data(9,2), kasset_c3_p25_data(9,2) !moment,wgt (inverse variance)
real(prec) :: passet_c3_p50_data(9,2), kasset_c3_p50_data(9,2) !moment,wgt (inverse variance)
real(prec) :: passet_c3_p75_data(9,2), kasset_c3_p75_data(9,2) !moment,wgt (inverse variance)
real(prec) :: passet_c4_p25_data(9,2), kasset_c4_p25_data(9,2) !moment,wgt (inverse variance)
real(prec) :: passet_c4_p50_data(9,2), kasset_c4_p50_data(9,2) !moment,wgt (inverse variance)
real(prec) :: passet_c4_p75_data(9,2), kasset_c4_p75_data(9,2) !moment,wgt (inverse variance)
real(prec) :: passet_c5_p25_data(7,2), kasset_c5_p25_data(7,2) !moment,wgt (inverse variance)
real(prec) :: passet_c5_p50_data(7,2), kasset_c5_p50_data(7,2) !moment,wgt (inverse variance)
real(prec) :: passet_c5_p75_data(7,2), kasset_c5_p75_data(7,2) !moment,wgt (inverse variance)
real(prec) :: FT_well_simu, FT_sick_simu
real(prec) :: informal_simu_1, informal_simu_2, informal_simu_avg 
real(prec) :: ins_simu_1, ins_simu_2, ins_simu_3, ins_simu_avg
real(prec) :: passet_w5_25_simu, passet_w5_50_simu, passet_w5_75_simu

integer :: simnode1,simnode2

!predata
real(prec), allocatable :: agesim96(:)
real(prec), allocatable :: HHIDsim(:)
real(prec), allocatable :: mstat96(:)
real(prec), allocatable :: asim96(:)
integer, parameter :: ncohort=5
integer, parameter :: nHRS=2433
integer, parameter :: nHRS1=327 
integer, parameter :: nHRS2=416 
integer, parameter :: nHRS3=651 
integer, parameter :: nHRS4=615 
integer, parameter :: nHRS5=424 
integer, parameter :: ntHRS=9 
integer HRSprof_age(nHRS), HRSprof_cohort(nHRS), HRSprof_ltci(nHRS)
double precision HRSprof_passets(nHRS), HRSprof_kassets(nHRS), HRSprof_pinc(nHRS), HRSprof_kinc(nHRS), HRSprof_health(nHRS,ntHRS)
integer HRSprof_age1(nHRS1), HRSprof_age2(nHRS2), HRSprof_age3(nHRS3), HRSprof_age4(nHRS4), HRSprof_age5(nHRS5) 
double precision HRSprof_passets1(nHRS1), HRSprof_passets2(nHRS2), HRSprof_passets3(nHRS3)
double precision HRSprof_passets4(nHRS4), HRSprof_passets5(nHRS5)
double precision HRSprof_kassets1(nHRS1), HRSprof_kassets2(nHRS2), HRSprof_kassets3(nHRS3)
double precision HRSprof_kassets4(nHRS4), HRSprof_kassets5(nHRS5)
integer HRSprof_pinc1(nHRS1), HRSprof_pinc2(nHRS2), HRSprof_pinc3(nHRS3), HRSprof_pinc4(nHRS4), HRSprof_pinc5(nHRS5)
integer HRSprof_kinc1(nHRS1), HRSprof_kinc2(nHRS2), HRSprof_kinc3(nHRS3), HRSprof_kinc4(nHRS4), HRSprof_kinc5(nHRS5)
integer HRSprof_health1(nHRS1,ntHRS), HRSprof_health2(nHRS2,ntHRS), HRSprof_health3(nHRS3,ntHRS)
integer HRSprof_health4(nHRS4,ntHRS), HRSprof_health5(nHRS5,ntHRS)



!real(prec),parameter,dimension(nm) :: paretoweight=(/0.0001,0.1,0.2,0.3,0.4,0.5,0.58,0.59,0.6,0.61,0.62,0.7,0.8,0.9,0.999/) !for joint pareto-altruism identification
real(prec),parameter,dimension(nm) :: paretoweight=(/0.0001,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,0.999/)
real(prec),parameter,dimension(nk) :: kgrid=       (/0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9/)

real(prec), parameter :: yrsinpd = 2.0 ! 1.0 !number of years in a period

! real(prec), parameter :: premium = 3170.3515*yrsinpd !0 load...
 real(prec), parameter :: premium = 3099.375*yrsinpd !yearly premium for LTC insurance (18% load)
! real(prec), parameter :: premium = 3630.6964*yrsinpd !yearly premium for LTC insurance (30% load)
! real(prec), parameter :: premium = 5082.975*yrsinpd !yearly premium for LTC insurance (50% load)
 real(prec), parameter :: premiumPO = 3099.375*yrsinpd !yearly premium for LTC insurance (18% load) PUBLIC OPTION

!hours spent
real(prec), parameter :: Tyear= 5840*yrsinpd !Number of available hours in a year (16*365)
real(prec), parameter :: workhrs(3)= (/ 0*yrsinpd,1000*yrsinpd,2000*yrsinpd /) !not working, part time, full time
real(prec), parameter :: famhrs(nh)= (/ 0*yrsinpd,1000*yrsinpd,2000*yrsinpd /) !no family care, full family care (based on health status)
real(prec), parameter :: meaninc(2)=(/ 24960.0d0*yrsinpd,45600.0d0*yrsinpd /) ! two kid income types
real(prec), parameter :: nonwageinc(nw)= 0d0

real(prec), parameter :: PTwagefrac=0.5

double precision hcost(ni,nh,2) !insurance state, health state, family care state
double precision mcdcash(nh,2) !health state, family care state

double precision disc, iniass
double precision reprate, minw, mw, rhow, r, minass

! Params for different processes and preferences
data minass	/ 1.d-10	/  ! (machine zero issue)
data disc	/ 1.064d0   /  !/ 1.0422 / !  (1+discount rate) -- 1.0422 is for annual periods

! Income Process
real(prec), parameter :: varw(2) = (/ 0.0277, 0.0437 /)
data rhow	    / 1.0d0    / ! Parameter measuring persistence in AR(1) process for income. For RW,=1. for iid,=0.
data mw         / 2.0d0    / ! multiple of standard deviations from mean for extremes of income(t)
data r			/ 1.03d0   / ! / 1.0149d0  / (1+interest rate) -- 1.0149 is for annual periods
data minw       / 0.0d0    / ! min income

double precision, allocatable :: survprob(:,:,:)
double precision, allocatable :: alowp(:), alowk(:), ahighp(:), ahighk(:)
double precision, allocatable :: assetp(:,:), assetk(:,:)
double precision, allocatable :: maxincp(:), maxinck(:), minincp(:), mininck(:) 
double precision, allocatable :: incleftp(:,:), incleftk(:,:)
double precision, allocatable :: hprob(:,:,:,:)

double precision, allocatable :: valueNT0_married(:,:,:,:,:,:),valueNT100_married(:,:,:,:,:,:)
double precision, allocatable :: valuePT0_married(:,:,:,:,:,:),valuePT100_married(:,:,:,:,:,:)
double precision, allocatable :: valueFT0_married(:,:,:,:,:,:),valueFT100_married(:,:,:,:,:,:)
double precision, allocatable :: value_married(:,:,:,:,:,:)

double precision, allocatable :: valuepNT0_married(:,:,:,:,:,:),valuepNT100_married(:,:,:,:,:,:)
double precision, allocatable :: valuepPT0_married(:,:,:,:,:,:),valuepPT100_married(:,:,:,:,:,:)
double precision, allocatable :: valuepFT0_married(:,:,:,:,:,:),valuepFT100_married(:,:,:,:,:,:)
double precision, allocatable :: valuep_married(:,:,:,:,:,:)

double precision, allocatable :: valuekNT0_married(:,:,:,:,:,:),valuekNT100_married(:,:,:,:,:,:)
double precision, allocatable :: valuekPT0_married(:,:,:,:,:,:),valuekPT100_married(:,:,:,:,:,:)
double precision, allocatable :: valuekFT0_married(:,:,:,:,:,:),valuekFT100_married(:,:,:,:,:,:)
double precision, allocatable :: valuek_married(:,:,:,:,:,:)

double precision, allocatable :: valuep_padj(:,:,:,:,:),valuek_padj(:,:,:,:,:)
double precision, allocatable :: valuep_kadj(:,:,:,:,:),valuek_kadj(:,:,:,:,:)

double precision, allocatable :: conexppNT0_married(:,:,:,:,:,:),conexppNT100_married(:,:,:,:,:,:)
double precision, allocatable :: conexppPT0_married(:,:,:,:,:,:),conexppPT100_married(:,:,:,:,:,:)
double precision, allocatable :: conexppFT0_married(:,:,:,:,:,:),conexppFT100_married(:,:,:,:,:,:)
double precision, allocatable :: conexpp_married(:,:,:,:,:,:)

double precision, allocatable :: conexpkNT0_married(:,:,:,:,:,:),conexpkNT100_married(:,:,:,:,:,:)
double precision, allocatable :: conexpkPT0_married(:,:,:,:,:,:),conexpkPT100_married(:,:,:,:,:,:)
double precision, allocatable :: conexpkFT0_married(:,:,:,:,:,:),conexpkFT100_married(:,:,:,:,:,:)
double precision, allocatable :: conexpk_married(:,:,:,:,:,:)

double precision, allocatable :: kappaNT0_married(:,:,:,:,:,:), kappaNT100_married(:,:,:,:,:,:)
double precision, allocatable :: kappaPT0_married(:,:,:,:,:,:), kappaPT100_married(:,:,:,:,:,:)
double precision, allocatable :: kappaFT0_married(:,:,:,:,:,:), kappaFT100_married(:,:,:,:,:,:)
double precision, allocatable :: kappa_married(:,:,:,:,:,:)


double precision, allocatable :: conexpp_padj(:,:,:,:,:),conexpk_padj(:,:,:,:,:)
double precision, allocatable :: conexpp_kadj(:,:,:,:,:),conexpk_kadj(:,:,:,:,:)


double precision, allocatable :: valuepNT_divorced(:,:,:,:,:),valuekNT_divorced(:,:,:,:,:)
double precision, allocatable :: valuepPT_divorced(:,:,:,:,:),valuekPT_divorced(:,:,:,:,:)
double precision, allocatable :: valuepFT_divorced(:,:,:,:,:),valuekFT_divorced(:,:,:,:,:)
double precision, allocatable :: valuep_divorced(:,:,:,:,:),valuek_divorced(:,:,:,:,:)

double precision, allocatable :: conexppNT_divorced(:,:,:,:,:),conexpkNT_divorced(:,:,:,:,:)
double precision, allocatable :: conexppPT_divorced(:,:,:,:,:),conexpkPT_divorced(:,:,:,:,:)
double precision, allocatable :: conexppFT_divorced(:,:,:,:,:),conexpkFT_divorced(:,:,:,:,:)
double precision, allocatable :: conexpp_divorced(:,:,:,:,:),conexpk_divorced(:,:,:,:,:)

double precision, allocatable :: LFP_divorced(:,:,:,:,:)
double precision, allocatable :: insurance_divorced(:,:,:)

real(prec), allocatable :: LFP_deadparent(:,:)
real(prec), allocatable :: valuekNT_deadparent(:,:), valuekPT_deadparent(:,:), valuekFT_deadparent(:,:), valuek_deadparent(:,:)
real(prec), allocatable :: conexpkNT_deadparent(:,:), conexpkPT_deadparent(:,:), conexpkFT_deadparent(:,:), conexpk_deadparent(:,:)

double precision, allocatable :: value(:,:,:,:,:,:),valuep(:,:,:,:,:,:),valuek(:,:,:,:,:,:)


double precision, allocatable :: consump(:,:,:,:,:,:),consumk(:,:,:,:,:,:)
double precision, allocatable :: kappa(:,:,:,:,:,:)
double precision, allocatable :: kappa_padj(:,:,:,:,:),kappa_kadj(:,:,:,:,:)
double precision, allocatable :: LFP_married(:,:,:,:,:,:)
double precision, allocatable :: LFP_kadj(:,:,:,:,:), LFP_padj(:,:,:,:,:)
double precision, allocatable :: LFP(:,:,:,:,:,:)
double precision, allocatable :: famcare_married(:,:,:,:,:,:)
double precision, allocatable :: famcare_kadj(:,:,:,:,:), famcare_padj(:,:,:,:,:)
double precision, allocatable :: family(:,:,:,:,:,:)
double precision, allocatable :: married(:,:,:,:,:,:)
double precision, allocatable :: mu_padj(:,:,:,:,:),mu_kadj(:,:,:,:,:)
double precision, allocatable :: mu(:,:,:,:,:,:)
double precision, allocatable :: insurance_married(:,:,:,:)
double precision, allocatable :: insurance_kadj(:,:,:),insurance_padj(:,:,:)
double precision, allocatable :: insurance(:,:,:,:)


real(prec), allocatable :: randincp(:),randinck(:),randincktype(:),randsurv(:),randhlth(:),randprofile(:) !random numbers
real(prec), allocatable :: randlovek(:),randlovep(:)

integer :: recsize1,recsize2,recsize3,recsize4,recsize5,recsize6,recsizesim,recsizesim2
real(prec), allocatable :: tempvec1(:), tempvec2(:), tempvec3(:),tempvec4(:), tempvec5(:), tempvec6(:),tempvecsim(:)
integer, allocatable :: tempvecsimint(:),tempvecsimint2(:)

contains
	!*****************************
	subroutine arrays	! allocate arrays
		
		implicit none
		
		nbignodes = nh*nwp
		nonodesk = nw*na  
		recsize1=nt*ni*nbignodes*nm*na*nonodesk
		recsize2=nt*ni*nbignodes*na*nonodesk
		recsize3=ni*nbignodes*nm*na*nonodesk
		recsize4=ni*nbignodes*na*nonodesk
		recsize5=nbignodes*nm*na*nonodesk
		recsize6=nbignodes*na*nonodesk
		recsizesim=nt*nosims
		recsizesim2=nosims
		
		allocate (tempvec1(recsize1),tempvec2(recsize2),tempvecsim(recsizesim),tempvec3(recsize3),tempvec4(recsize4), &
				     tempvec5(recsize5),tempvec6(recsize6),&
			      tempvecsimint(recsizesim),tempvecsimint2(recsizesim2),&
					agesim96(nosims),HHIDsim(nosims),mstat96(nosims),asim96(nosims),&
						survprob(nh,nwp,nt), &
					ahighp(nt),ahighk(nt),alowp(nt),alowk(nt), &
				  incleftp(nwp,nt), incleftk(nw,nt), &
				  maxincp(nt), maxinck(nt), minincp(nt), mininck(nt), &
				  assetp(na,nt), assetk(na,nt), &
				  hprob(nh,nh,nwp,nt), &

				  valuekNT_deadparent(nt,nonodesk),valuekPT_deadparent(nt,nonodesk),valuekFT_deadparent(nt,nonodesk), &
				  valuek_deadparent(nt,nonodesk), &
				  conexpkNT_deadparent(nt,nonodesk),conexpkPT_deadparent(nt,nonodesk),conexpkFT_deadparent(nt,nonodesk), &
				  conexpk_deadparent(nt,nonodesk), &	
				  LFP_deadparent(nt,nonodesk), &			  
				  
				  valuepNT0_married(nt,ni,nbignodes,nm,na,nonodesk),valuekNT0_married(nt,ni,nbignodes,nm,na,nonodesk), &
				  valuepNT100_married(nt,ni,nbignodes,nm,na,nonodesk),valuekNT100_married(nt,ni,nbignodes,nm,na,nonodesk), &
				  valuepPT0_married(nt,ni,nbignodes,nm,na,nonodesk),valuekPT0_married(nt,ni,nbignodes,nm,na,nonodesk), &
				  valuepPT100_married(nt,ni,nbignodes,nm,na,nonodesk),valuekPT100_married(nt,ni,nbignodes,nm,na,nonodesk), &
				  valuepFT0_married(nt,ni,nbignodes,nm,na,nonodesk),valuekFT0_married(nt,ni,nbignodes,nm,na,nonodesk), &
				  valuepFT100_married(nt,ni,nbignodes,nm,na,nonodesk),valuekFT100_married(nt,ni,nbignodes,nm,na,nonodesk), &
				  valuep_married(nt,ni,nbignodes,nm,na,nonodesk),valuek_married(nt,ni,nbignodes,nm,na,nonodesk), &

					valueNT0_married(nt,ni,nbignodes,nm,na,nonodesk), valueNT100_married(nt,ni,nbignodes,nm,na,nonodesk), &
					valuePT0_married(nt,ni,nbignodes,nm,na,nonodesk), valuePT100_married(nt,ni,nbignodes,nm,na,nonodesk), &
					valueFT0_married(nt,ni,nbignodes,nm,na,nonodesk), valueFT100_married(nt,ni,nbignodes,nm,na,nonodesk), &
					value_married(nt,ni,nbignodes,nm,na,nonodesk), &

				  conexppNT0_married(nt,ni,nbignodes,nm,na,nonodesk),conexpkNT0_married(nt,ni,nbignodes,nm,na,nonodesk), &
				  conexppNT100_married(nt,ni,nbignodes,nm,na,nonodesk),conexpkNT100_married(nt,ni,nbignodes,nm,na,nonodesk), &
				  conexppPT0_married(nt,ni,nbignodes,nm,na,nonodesk),conexpkPT0_married(nt,ni,nbignodes,nm,na,nonodesk), &
				  conexppPT100_married(nt,ni,nbignodes,nm,na,nonodesk),conexpkPT100_married(nt,ni,nbignodes,nm,na,nonodesk), &
				  conexppFT0_married(nt,ni,nbignodes,nm,na,nonodesk),conexpkFT0_married(nt,ni,nbignodes,nm,na,nonodesk), &
				  conexppFT100_married(nt,ni,nbignodes,nm,na,nonodesk),conexpkFT100_married(nt,ni,nbignodes,nm,na,nonodesk), &
				  conexpp_married(nt,ni,nbignodes,nm,na,nonodesk),conexpk_married(nt,ni,nbignodes,nm,na,nonodesk), &

					kappaNT0_married(nt,ni,nbignodes,nm,na,nonodesk), kappaNT100_married(nt,ni,nbignodes,nm,na,nonodesk), &
					kappaPT0_married(nt,ni,nbignodes,nm,na,nonodesk), kappaPT100_married(nt,ni,nbignodes,nm,na,nonodesk), &
					kappaFT0_married(nt,ni,nbignodes,nm,na,nonodesk), kappaFT100_married(nt,ni,nbignodes,nm,na,nonodesk), &
					kappa_married(nt,ni,nbignodes,nm,na,nonodesk), &

				  valuep_padj(nt,ni,nbignodes,na,nonodesk),valuek_padj(nt,ni,nbignodes,na,nonodesk), &
				  valuep_kadj(nt,ni,nbignodes,na,nonodesk),valuek_kadj(nt,ni,nbignodes,na,nonodesk), &
				  valuepNT_divorced(nt,ni,nbignodes,na,nonodesk),valuekNT_divorced(nt,ni,nbignodes,na,nonodesk), &
				  valuepPT_divorced(nt,ni,nbignodes,na,nonodesk),valuekPT_divorced(nt,ni,nbignodes,na,nonodesk), &
				  valuepFT_divorced(nt,ni,nbignodes,na,nonodesk),valuekFT_divorced(nt,ni,nbignodes,na,nonodesk), &
				  valuep_divorced(nt,ni,nbignodes,na,nonodesk),valuek_divorced(nt,ni,nbignodes,na,nonodesk), &
				  conexppNT_divorced(nt,ni,nbignodes,na,nonodesk),conexpkNT_divorced(nt,ni,nbignodes,na,nonodesk), &
				  conexppPT_divorced(nt,ni,nbignodes,na,nonodesk),conexpkPT_divorced(nt,ni,nbignodes,na,nonodesk), &
				  conexppFT_divorced(nt,ni,nbignodes,na,nonodesk),conexpkFT_divorced(nt,ni,nbignodes,na,nonodesk), &
				  conexpp_divorced(nt,ni,nbignodes,na,nonodesk),conexpk_divorced(nt,ni,nbignodes,na,nonodesk), &
				  value(nt,ni,nbignodes,nm,na,nonodesk), &
				  valuep(nt,ni,nbignodes,nm,na,nonodesk),valuek(nt,ni,nbignodes,nm,na,nonodesk), &				  
				  conexpp_padj(nt,ni,nbignodes,na,nonodesk),conexpk_padj(nt,ni,nbignodes,na,nonodesk), &
				  conexpp_kadj(nt,ni,nbignodes,na,nonodesk),conexpk_kadj(nt,ni,nbignodes,na,nonodesk), &
				  consump(nt,ni,nbignodes,nm,na,nonodesk),consumk(nt,ni,nbignodes,nm,na,nonodesk), &
				  married(nt,ni,nbignodes,nm,na,nonodesk), &
				  kappa(nt,ni,nbignodes,nm,na,nonodesk), &
				  kappa_padj(nt,ni,nbignodes,na,nonodesk), kappa_kadj(nt,ni,nbignodes,na,nonodesk), &
				  LFP_divorced(nt,ni,nbignodes,na,nonodesk), &
				  LFP_married(nt,ni,nbignodes,nm,na,nonodesk), &
				  LFP_kadj(nt,ni,nbignodes,na,nonodesk), LFP_padj(nt,ni,nbignodes,na,nonodesk), &
				  LFP(nt,ni,nbignodes,nm,na,nonodesk), &
				  famcare_married(nt,ni,nbignodes,nm,na,nonodesk), &
				  famcare_kadj(nt,ni,nbignodes,na,nonodesk), famcare_padj(nt,ni,nbignodes,na,nonodesk), &
				  family(nt,ni,nbignodes,nm,na,nonodesk), &
				  mu_padj(nt,ni,nbignodes,na,nonodesk), mu_kadj(nt,ni,nbignodes,na,nonodesk), &
				  mu(nt,ni,nbignodes,nm,na,nonodesk), &
				  insurance_divorced(nbignodes,na,nonodesk), &
				  insurance_married(nbignodes,nm,na,nonodesk), &
				  insurance_kadj(nbignodes,na,nonodesk), insurance_padj(nbignodes,na,nonodesk), &
				  insurance(nbignodes,nm,na,nonodesk), &

				  randincp(nosims*(nt)),randinck(nosims*(nt)),randincktype(nosims),randlovep(nosims*(nt)),randlovek(nosims*(nt)),&
						randsurv(nosims*(nt)),randhlth(nosims*(nt)),randprofile(nosims*(nt))&
				 )

	!health cost for healthy, semi-unhealthy, unhealthy: insurance state, health state, family care state
	!hcost(ni,nh,2) !insurance state, health state, family care state
	if (simcashins==1) then 
		hcost(1,:,1)=(/ 0.0D0,-20000.0D0*dble(yrsinpd),-44350.0D0*dble(yrsinpd) /) !has insurance (61700=cost, insurance pays max 44350 so still have to pay OOP 17350)
		hcost(1,:,2)=(/ 0.0D0,     0.0D0*dble(yrsinpd), 17350.0D0*dble(yrsinpd) /) !has insurance (61700=cost, insurance pays max 44350 so still have to pay OOP 17350)
		hcost(2,:,1)=(/ 0.0D0,     0.0D0*dble(yrsinpd),     0.0D0*dble(yrsinpd) /) !does not have insurance
		hcost(2,:,2)=(/ 0.0D0, 20000.0D0*dble(yrsinpd), 61700.0D0*dble(yrsinpd) /) !does not have insurance
	elseif (simfullfairins==1) then 
		hcost(1,:,1)=(/ 0.0D0,     0.0D0*dble(yrsinpd),     0.0D0*dble(yrsinpd) /) !has insurance (61700=cost, insurance pays max 44350 so still have to pay OOP 17350)
		hcost(1,:,2)=(/ 0.0D0,     0.0D0*dble(yrsinpd),     0.0D0*dble(yrsinpd) /) !has insurance (61700=cost, insurance pays max 44350 so still have to pay OOP 17350)
		hcost(2,:,1)=(/ 0.0D0,     0.0D0*dble(yrsinpd),     0.0D0*dble(yrsinpd) /) !does not have insurance
		hcost(2,:,2)=(/ 0.0D0, 20000.0D0*dble(yrsinpd), 61700.0D0*dble(yrsinpd) /) !does not have insurance
	elseif (simpublicoption==1) then
		hcost(1,:,1)=(/ 0.0D0,     0.0D0*dble(yrsinpd),     0.0D0*dble(yrsinpd) /) !has insurance (61700=cost, insurance pays max 44350 so still have to pay OOP 17350)
		hcost(1,:,2)=(/ 0.0D0,     0.0D0*dble(yrsinpd), 17350.0D0*dble(yrsinpd) /) !has insurance (61700=cost, insurance pays max 44350 so still have to pay OOP 17350)
		hcost(2,:,1)=(/ 0.0D0,     0.0D0*dble(yrsinpd),     0.0D0*dble(yrsinpd) /) !does not have insurance
		hcost(2,:,2)=(/ 0.0D0, 20000.0D0*dble(yrsinpd), 61700.0D0*dble(yrsinpd) /) !does not have insurance				
		hcost(3,:,1)=(/ 0.0D0,-20000.0D0*dble(yrsinpd),-44350.0D0*dble(yrsinpd) /) !has insurance public option (61700=cost, insurance pays max 44350 so still have to pay OOP 17350)
		hcost(3,:,2)=(/ 0.0D0,     0.0D0*dble(yrsinpd), 17350.0D0*dble(yrsinpd) /) !has insurance public option (61700=cost, insurance pays max 44350 so still have to pay OOP 17350)	
	else 	
		hcost(1,:,1)=(/ 0.0D0,     0.0D0*dble(yrsinpd),     0.0D0*dble(yrsinpd) /) !has insurance (61700=cost, insurance pays max 44350 so still have to pay OOP 17350)
		hcost(1,:,2)=(/ 0.0D0,     0.0D0*dble(yrsinpd), 17350.0D0*dble(yrsinpd) /) !has insurance (61700=cost, insurance pays max 44350 so still have to pay OOP 17350)
		hcost(2,:,1)=(/ 0.0D0,     0.0D0*dble(yrsinpd),     0.0D0*dble(yrsinpd) /) !does not have insurance
		hcost(2,:,2)=(/ 0.0D0, 20000.0D0*dble(yrsinpd), 61700.0D0*dble(yrsinpd) /) !does not have insurance			
	endif

	!mcdcash(nh,2) !health state, family care state
	mcdcash=0.0D0
	if (simcashmcd==1) then
		mcdcash(:,1)=(/ 0.0D0,20000.0D0*dble(yrsinpd),61700.0D0*dble(yrsinpd) /) !using family care
	endif

	end subroutine arrays

end module params
